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ABSTRACT 

We run very large cosmological iV-body hydrodynamical simulations in order to 
study statistically the baryon fractions in early dark matter halos. We critically exam- 
ine how differences in the initial conditions affect the gas fraction in the redshift range 
z = 11 — 21. We test three different linear power spectra for the initial conditions: (1) A 
complete heating model, which is our fiducial model; this mod el follows the evolution 
of overdensities correctly, according to lNaoz fc Barkana] ([2005), in particular including 
the spatial variation of the speed of sound of the gas due to Compton heating from 
the CMB. (2) An equal-5 model, which assumes that the initial baryon fluctuations 
are equal to those of the dark matter, while conserving ag of the total matter. (3) A 
mean c s model, which assumes a uniform speed of sound of the gas. The latter two 
models are often used in the literature. We calculate the baryon fractions for a large 
sample of halos in our simulations. Our fiducial model implies that before reionization 
and significant stellar heating took place, the minimum mass needed for a minihalo 
to keep most of its baryons throughout its formation was ~ 3 x 10 4 M . However, 
the alternative models yield a wrong (higher by about 50%) minimum mass, since the 
system retains a memory of the initial conditions. We also demonstrate this using the 
"filtering mass" from linear theory, which accurately describes the evolution of the 
baryon fraction throughout the simulated redshift range. 



• • 1 INTRODUCTION 



Recent measurements of anisotropies of the cosmic mi- 
crowave background (CMB) radiation have revealed the de- 
tailed distribution of matter in th e Universe a few hundred 
thousand years after t h e Big Bang |Spergel et al 1l200il2007l; 
iKomatsu etafl 120081 . I2010T ). Observations utilizing large 
ground-based telescopes and space telescopes have discov- 
ered galaxies and black holes that were in place when the age 
of the Universe was less than a billion years. Moreover, many 
galaxies have been found at z > 7 (|Bouwens et alj l20ld : 
iMcLure et aliioioh in the Hubble Ultra Deep Field, whereas 
already a few gamma-r ay bursts at z > 6 h a ve been detected 
by the Swift satellite jTanyir et all 120091 ; ISalvaterra et al.l 
120091 ; iLin. Liang, fc Zhand l201Ch . These first objects are 
probably the building blocks of the present day galaxies, 
thus, solving the puzzle behind their formation will have a 
profound implication o n our understanding of the Universe 
(see for recent reviews iBromm et all |2009| ; lYoshidal 120091, 
and references therein). 

The formation of the first generation of galaxies in the 
Universe has been studied for many years. High resolution 
cosmological simulations can follow complex astrophysical 
processes, while analytical calculations can provide an over- 
all understanding, and can be used to decouple different 



physical effects seen in simulations. Analytic models are also 
useful for estimating the limitations of numerical simulations 
such as insufficient resolution and small boxsizes (Yoshida 
et al. 2003; Barkana & Loeb 2004; Naoz & Barkana 2005). 
Combining the two approaches may offer many of the ad- 
vantages of both. 

The initial conditions (hereafter ICs) in a cosmologi- 
cal simulation can have a large effect on the formation of 
the first galaxies in simulations, i.e., both on the formation 
time (or on the halo abundance at a given time) and the 
halo prop erties at formation t ime (such as the average gas 
fraction). lYoshida et all (|2003n studied high-redshift struc- 
ture formation and reionization while testing two different 
models for power spectra as their ICs. They found that 
different models have a profound effect on the abundance 
of primordial star-forming gas clouds and thus on when 
the reionization was i nitiated and its progress . In th e an- 
alytical point of view. iNaoz. Noter &: Barkanal (|2006T ) and 
iNaoz fc Barkanal l|2007l ) showed that the ICs at high red- 
shift have a significant effect on the halo abundance and the 
gas fraction at virialization. While these effects are largest 
at the hi ghest redshift, e.g., z ~ 65 for the first star in the 



universe ( Naoz. Noter fc Bark ana 2006), they are still sig- 
nificant for halos forming at z ~ 10 — 30. The first gas rich 
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halos at these redshifts are expected to hos t the first stars 
(z ~ 65 - 30 iNaoz. Noter fc Barkanal 120061; lYoshidal l200rj ; 
iGao et al l 120071 : iTrenti. Stiavelli. fe Michael Shulll l2009al) 
and even the first gamma-ray bursts (e.g.. iBromm fc Loebl 
l2006l ; lNaoz fe Brombergil2007l ). Thus, investigating the for- 
mation properties of these halos is of prime importance . 

Gas rich halos in the early Universe may very well 
be a nurturing ground for dwar f galaxies, which at high 
redshift can form stars (e.g., Bromm. C oppi. & Larson 



2002, 19 991; lAbel. Bryan, fc Normanl 120021; lYoshida et al 



20061 ; lYoshida. Omukai fc HernauistJ 
ences therein) perhaps even at a hi 
rate 



200. 



and refer- 
star formation 



tncrcm) pcrriaps even a t a ni gh star tormation 
(iRicotti, Gnedin, fc Shuil |2002| ; iGreif et all |2010| ; 
IClark et al.l l2010h . Their properties are very important as 
they are responsible for metal poll ution and the ionizing 
radia t ion at these early t imes (e.g.. IShapiro. Iliev. fc Ragal 
20041; ICiardi et al.l 120061; iGnedin. Kravtsov fc Chenl |200S| ; 
Trenti fc Stiavellil |2009b). Moreover, halos that are too 
small for efficient cooling via atomic hydrogen, i.e., 
minihalos, are most susceptible to the effect of initial 
conditions. While they may not normally host astrophys- 
ical sources, minihalos may pro d uce a 21 -cm signature 



Kuhlen. Madau. fc Montgomery! (|2006t ) ; IShapiro et al 
200e ); INaoz fc Barkanal IMPOST but see iFurlanetto fc Oh 
2006)), and they can block ionizing radiation and 



produce an overall delay in the initial progress of reion- 
ization (e.g.. iBarkana fc Loebl bool ; llliev et al.ll2003l . 120051 ; 
iMcQuinn et al.ll2007h . The evolution of the halo gas fraction 
at various epochs of the universe is of prime importance, 
particularly in the early universe. 

In this paper, we examine the effect of using different 
initial conditions in simulations on the resulting minimum 
gas-rich halo ma ss in the redshift regime z = 11 — 21 . We per- 
form Gadget-2 (|Springel. Yoshida. fc White] l200ll ; ISpringell 
120051 ) simulations using a total of 768 3 x 2 particles. We 
compa re the initial conditions presented in lNaoz fc Barkanal 
(2005), which describe the linear evolution of overdensities 
in a fully consistent way, to two other alternative ICs, often 
used in the literature. We also compare to the prediction of 
the gas-rich mass from linear theory. We describe our dif- 
ferent initial conditions and simulations in sections [2] and 
[3] respectively. Our simulation results are presented in sec- 
tion U where we divide our discussion to the evolution of 
the non- linear power spectra (section |4.1[) and to the mini- 
mum gas-rich halo mass resulting from either linear theory 
or from the simulations (section l4.2p . Finally, we discuss our 
conclusions (section [5]). 

Throughout this paper, we adopt the following cosmo- 



logical parameters: (O.A, J2m, f2f>, 
0.046, 1, 0.82, 70 km s _1 Mpc" 1 ) 



n, q 8 , H )= (0 72, .28, 
|Komatsu et al.ll2009l ). 



2 DIFFERENT INITIAL CONDITION 
MODELS - BASIC EQUATIONS 

2.1 The fiducial ICs - "fid" 



We follow INaoz fc Barkanal l|2005h . who studied the linear 
evolution of both dark matter and baryon overdensities. The 
fluctuations of the temperature of the baryons (St) can- 
not be described as a simple function of a spatially uniform 
baryonic sound speed c 3 (t), as was previously assumed (e.g., 
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Figure 1. The relative difference (specifically, S mod ^i/8 ch — 1) 
between the fiducial linear initial conditions and the alternative 
models at z = 99. We consider the relative difference between the 
fid ICs and the mean c s ICs for both the baryons and dark matter 
(solid and short-dashed curves, respectively), and the relative dif- 
ference between the fid ICs and the E-5 ICs for both the baryons 
and dark matter (dotted and long-dashed curves, respectively). 
Note that we have plotted here the absolute value; the mean c s 
model gave a negative value (i.e., an underestimate compared to 
the fid model) while the E-<5 model gave a positive value (i.e., an 
overestimate) . 



iMa fc Bertschingerlll995l ). Furthermore, at high redshifts, 
the baryon density fluctuations (<5b) are not equal to those 
of dark matter (5dm) (contrary to a common assumption in 
simulations; four redsh ift examples are shown in figure 1 of 
INaoz fc Barkanall2005h . We label the power spectrum model 
as the "fid" (fiducial) ICs since it follows the evolution of lin- 
ear overdensiti es in a complete and con sistent way. 

Following INaoz fc Barkanal (|2005h we write the basic 
equations that describe the evolution of the dark matter, 
baryon density and temperature fluctuations: 



<5dm + 2ff<5dm — TJ-ffo — J- (fbSb + /dm<5dm) 



(1) 



where /dm and /b are the mean cosmic dark matter and 
baryonic fraction respectively. Here we follow the standard 
notations for cosmological parameters such as f2 m , Ho. The 
baryons are also subject to a pressure term: 



5 b +2H5 b = ^H^ 
2 



(fbSb + fdmSd 



k 2 k B T 
a 2 fi 



(5b + S T ) , (2) 



where jj, is the mean molecular weight, fcg is the Boltzmann 
constant and k i s the wavenumbe r . Usin g the first law of 
thermodynamics, INaoz fc Barkanal (|2005h derived the equa- 
tions for the evolution of the baryon average temperature 
and temperature fluctuations: 



dT 
It 



T) a 



(3) 
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where T 7 = [2.725 K]/a is the mean CMB temperature, and 
the first-order equation for the perturbation: 

dS T 2d8 b x e (t) _ 4 f r (% X %. . xl lA . 

with the second term on the right-hand side accounting for 
the Compton scattering of the CMB photons on the residual 
electrons from recombination, where x e (t) is the electron 
fraction out of the total number density of gas particles at 
time t, and 

J.-1 — 8_ (TtC m-13 -1 /r\ 

ty = -p~f = 8.55 x 10 yr , (5) 

where <jt is the Thomson scattering cross section and p 1 is 
the photon energy density. The first term on the right-hand- 
side of each of these two equations ((3]) and ([4]) accounts for 
adiabatic expansion of the gas, and the remaining terms 
capture t he effect of the thermal exchange with the CMB. 
Following IN aoz fc Barkanal ([2005) we have numerically cal- 
culated the evolut ion of the perturbations by modifying the 
CMBFAST code (jSeliak fc Zaldarriaeal 119961 ) according to 
th ese equations. Note that simila r physics was also explored 
bv lYamamoto et al l (|l997l . ll998l ). 

We solve the complete set of equations to obtain the 
power spectrum at different redshifts which can be used as 
initial conditions for our simulations. Figure [T] shows the 
ratio between this initial condition to the two alternative 
models tested in this paper. 



2.2 Alternative model I - equal 8 - "E-<5" 

In many cosmological ICs for N-body simulations and semi- 
analytical calculations, the fluctuations of the baryons are 
assumed to be equal to the fluctuations of the dark matter. 
We construct a model that includes this incorrect assump- 
tion while maintaining the correct overall Stot (i.e., conserv- 
ing as at z — 0, see appendix 1X1 for more details). Thus, in 
our "E-<5" model we calculate the correct Stot as a combina- 
tion of Sh and o~d m from the fiducial calculation in section 
(|2.1|l , but then take the baryon perturbation to be the same 
as for the dark matter, namely: 

<rE5 zES cch e <-ch . r srch ff*\ 
Ob = Odm = O to t = 7bO b + /dm Odm , \P) 

where 8^ h dul (o"jf dm ) is the resulting linear over-density from 
the fiducial calculation (E-<5 model) for the baryons and dark 
matter, respectively. We then compare the equal 8 model to 
our fiducial calculation. Figure [T] shows the ratio between 
the fid ICs and the E-<5 model for both the baryons and dark 
matter. We find that the E-<5 model overestimates the baryon 
fluctuations by > 30% on large scales (fc -1 > 10 kpc) while 
the overestimate grows to a much larger factor on small 
scales. 

Before recombination the baryons were tightly coupled 
to the radiation, resulting in suppression of the growth of 
their overdensity. However, the dark matter component, 
which is not affected by the photons, could basically grow 
once the fluctuation wavelength entered the Hubble horizon 
(in the linear regime, before equality, the dark mater fluctu- 
ations grew logarithmically with the scale factor, where after 
equality they grew linearly with the scale factor). Therefore, 
this resulted in a suppression of the baryonic overdensity by 



about three orders of magnitude co mpare to the dark mat - 
ter at recombination (e.g., fig. 1 in iNaoz fc Barkanal lioblih . 
While the baryons subsequently fall into the potential wells 
of the dark matter, it takes them some time to catch up, and 
the baryon fluctuations are still suppressed even at lower 
redshifts. This point is often overlooked in simulations and 
analytical calculations. 

2.3 Alternative model II - the mean sound speed 
approximation - "mean cs" 

INaoz fc Barkanal l|2005h showed that the presence of spa- 
tial fluctuations in the sound speed modifies the calculation 
of perturbation growth significantly. Nevertheless, for com- 
pleteness and as a case of comparison with previous results, 
we compare the simulation results with the results obtained 
using this approximation. Thus, we proceed by presenting 
the basic equations of the growth of density fluctuations, 
in this approximation of a uniform sound speed (hereafter 
"mean c s "). The evolution of the density fluctuations is de- 
scribed by a different set of coupled second order differential 
equations: 

3 2 ^2 

8dm + 2H8dm = 5- (/bOb + /dm Odm) , (7) 

2 a d 

Ob + 2H5 h = ^rHl — ^ (/bOb + /dm 5dm) — ^-CsOb , 

where c a = dp/dp is assumed to be spatially uniform (i.e., 
independent of k) and is thus calculated from the thermal 
evolution of a uniform gas undergoing Hubble expansion. 
With this assumption, the temperature fluctuations (as a 
function of k) are simply proportional at any given time to 
the gas density fluctuations: 



5b ksT I p 

INaoz fe Barkanal (|2005h showed that this approximation 
leads to an underestimation of the baryon density fluctu- 
ations by up to 30% at z = 100 and 10% at z = 20 for large 
wavenumbers. Figure [1] shows the ratio between the mean c a 
initial conditions and the fiducial ones for both the baryons 
and dark matter. It agrees with our previous results, showing 
that the underestimate by the mean c s model is greatest at 
k' 1 ~ 1 kpc. The non-linear evolution resulting from these 
initial conditions will result in shallower potential wells com- 
pared to the fiducial calculation, 

Even though it is clear that the precise baryon temper- 
ature fluctuations at high redshift are very significant, still 
many simulations use initial conditions that assume a uni- 
form speed of sound in the Universe. As shown below this 
leads to significantly different estimates for the gas content 
of the early halos. 



3 THE SIMULATION 
3.1 Basic parameters 

We r u n a Gadget 2 s imulation (|Springel. Yoshida. fc White] 
l200ll : ISpringel Booll ) starting from redshift 99, for a to- 
tal of 2 x 768 3 particles (768 3 particles each for the Dark 
Matter and baryon components) and our box size is: 
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2 Mpc. We choose this box size so that a halo mass of 
10 5 Ma would have ~ 500 p articl es. This way according to 
iNaoz. Barkana fc Mesingerl ((2009) we are able to estimate 
the gas fraction in ~ 10 5 Mq halos correctly (see below for 
the halo definition). Our softening length is 0.2 comoving 
kpc. 

For all runs, glass-like cosmological ICs were gen- 
erated using the Zel'dovich approximation. The transfer 
functions were generated using the various models de- 
scribed above. We have used a glass file which was ran- 
domly displaced thus removing the coupling between nearby 
DM and gas particles. Using this randomization procedure 
we achieve essentially the same effect to that shown in 
lYoshida. Sugivama. fc Hernquistl (|2003l ). In generating the 
ICs, a convolution between the glass file and the transfer 
function from the different models was done, thus taking 
into account the different velocities of the DM and baryons 
(for the fiducial and mean cs models). We note that we have 
used the same phases for the DM and baryons, in all of the 
simulations. 

We set the initial temperature to be 164. 1 1 K (as de- 
rived from linear theory), and thus Gadget assumes neutral 
and monoatomic gas, and converts to thermal energy (i.e., 
adiabatic initial conditions) . Although this work emphasizes 
the need for a precise calculation of the baryon overdensities 
resulting from temperature fluctuations, we actually neglect 
the temperature fluctuations in the initial conditions. This 
may not be a bad approximation since the halos we study 
are already somewhat non-linear at our initial redshift, and 
the Compton heating is quite small compared to the adi- 
abatic heating during non-linear gravitational collapse (see 
Appendix [Cj. A more complete treatment would be to in- 
clude in the simulation the precise temperature fluctuations, 
which we leave for future work. Nevertheless, even with the 
current treatment our results show consistency with linear 
theory. 



3.2 Halo definition 

We locate dark matter halos by running a FOF group-finder 
algorithm with a linking parameter of 0.2. We then find the 
center of mass of each halo and calculate the density profile 
of the dark matter and baryons, separately. In order to derive 
the density profile we assume a spherical halo, and divide it 
to 2000 shells. Combining these density profiles, we find the 
virial radius rw at which the overdensity is 200 times the 
background density, and the gas fraction of each halo. 

Recently, iTrenti et alj j20f OS ) performed a resolution 
analysis in order to study the mass definition of halos in 
simulations. Their conclusion (their figure 2) is that using 
the FOF algorithm and assuming about 500 particles per 
spherical halo introduces an error of ~ 15% in the mass 
definition. In our gas fraction analysis we have chosen only 
halos with a number of particles larger or equal to 500, i.e., 
we limit our errors in halo mass definition to b e low ~ 15%. 
Also, according to INaoz. Barkana fc Mesingerl ((2009) , this 
way we can estimate the gas fraction inside a halo accu- 
rately. 
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Figure 2. Comparison of the linear and non-linear power 
spectra. The linear power spectra (generated according to 
INaoz fc Barkana 2005) are shown for the dark matter and baryon 
components (solid and dashed curves, respectively), while the cor- 
responding non-linear spectra (as measured in the simulations) 
are shown as triangles and squares, respectively. We show results 
for the fid model at redshifts z = 99, 21 and 11. 

4 RESULTS AND COMPARISON AMONG 
THE MODELS 

4.1 Non-linear power spectrum evolution 

One way to probe cosmic structure particularly on small 
scales is through the non-linear power spectrum. We be- 
gin our simulation at z = 99 with linear initial conditional 
The main disagreement between the three models lies in the 
baryonic component (although the E-<5 calculation also un- 
derestimates the dark matter overdensities by ~ 10% ). This 
input differenc e is then modified by the non-linea r evolu tion. 

Following lYoshida. Sugivama. fc Hernauistl (|2003T ) we 
compared the linear power spec t rum for the fid model, as 
computed from lNaoz fc Barkanal {2005), for the dark matter 
and baryon components, with the non-linear power spectra 
from the simulation (see figure [2}. The two power spectra 
agree well as expected in the linear regime. We note that 
the other two models approach the fid model at low red- 
shifts (see appendix IA1 figure I Al[l . 

FigElshows the differences among the fid, E-<5 and mean 
c s ICs, in terms of the non-linear power spectra at the later 
redshifts at which halos were formed in our simulation. The 
mean c s model maintains over time roughly the same level of 
discrepancy with the fid model, while in the E-<5 model both 
the baryonic and dark matter differences decline slightly 
slower than with the inverse scale factor. As clearly can 

1 This is, of course, a n app roximation, , since as shown in 
INaoz. Noter fc Barkanal d2006T ) at z = 99 overdensities are al- 
ready slightly non-linear. The effect of starting the simulation at 
high redshifts is studied elsewhere (Naoz et. al, in prep.) 
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Figure 3. The ratio of the non-linear power spectra (specifically, 
^model/^fid — 1) at 2 = 21, 15, and 12 (from bottom to top); 
curves are denoted as in figure [T] Note that we have plotted here 
the absolute value; the mean c s model underestimates and the 
E-5 model overestimates the power spectrum compared to the fid 
model. 



be seen from figure [3] the non-linear evolution of halos is 
still strongly affected by th e choice of initial condit ions even 
at redshift 12. The fid ICs l|Naoz fc Barkanall2005h describe 
the linear evolution consistently and thus represent the best 
available prescription for the initial conditions. 



4.2 The minimum gas rich mass 

Studying the galaxy evolution and reionization either by us- 
ing simulations (both AMR and SPH) or by using analytical 
calculations relies on knowing the amount of gas within the 
dark matter halos. The simplest assumption, often used in 
the literature, is that a dark matter halo has the mean cos- 
mic fraction. This can lead to incorrect results, especially 
when one tries to study star formation, galaxy mergers, and 
related phenomena. 

Consider the various scales involved in the formation of 
non-linear objects containing DM and gas. On large scales 
(small wavenumbers) gravity dominates halo formation and 
gas pressure can be neglected. On small scales, on the other 
hand, the pressure dominates gravity and prevents baryon 
density fluctuations from growing together with the dark 
matter fluctuations. The relative force balan ce at a given 
time can be characterized by the I Jeansl l|l92Sl ) scale, which 
is the minimum scale on which a small gas perturbation will 
grow due to gravity overcoming the pressure gradient. As 
long as the Compton scattering of the CMB on the residual 
free electrons after cosmic recombination kept the gas tem- 
perature coupled to that of the CMB, the Jeans mass was 
constant in time. However, at z ~ 200 the gas temperature 
decoupled from the CMB temperature and the Jeans mass 
began to decrease with time as the gas cooled adiabatically. 



Any overdensity on a scale more massive than the Jeans 
mass at a given time can begin to collapse, due to a lack of 
sufficient pressure. However, the Jeans mass is related only 
to the evolution of perturbations at a given time. When the 
Jeans mass itself varies with time, the overall suppression 
of the growth of perturbations depends on a time-averaged 

Jea ns mass. 

iGnedin fc Huil (19 981 defined a "filtering mass" that de- 
scribes the highest mass scale on which the baryonic pressure 
still mana ges to suppress the linear baryonic fluctuations sig- 
nificantly. IGnedin! (2000) suggested, based on a simulation, 
that the filtering mass also describes the largest halo mass 
whose gas content is significantly suppressed compared to 
the cosmic baryon fraction. The latter mass scale, in gen- 
eral termed the "characteristic mass" , is defined as the halo 
mass for which the enclosed baryon fraction equals half the 
mean cosmic fraction. Thus, the characteristic mass distin- 
guishes between gas-rich and gas-poor halos. Many semi- 
analytical models of dwarfs galaxies often use the character- 
istic ma ss scale in order to estimate the gas fraction in ha- 
los (e.g ., iBullock. Kravtsov. fc Weinberdl2000l ; H enson et al.l 
2002a, b; Somerville 2002). Theoretically this sets an approx- 
imate minimum value on the mass that can still form stars. 



4-2.1 Prediction from linear theory 

In linear the o ry th e filtering mass, first defined by 
IGnedin fc Hul 0998), describes the highest mass scale 
on which the baryon density fluctuations are sup- 
pressed s ignificantly comp a red t o the dark matter fluc- 
tuations. iNaoz fc Barkanal ((2007) included the fact that 
the baryons have smo other ICs than the dark matter (see 
INaoz fc Ba rkana 200§j) and found a lower value of the filter- 
ing mass (by a factor of 3 — 10, d epending on the redshift). 
Following INaoz fc Barkanal (120071 ). the filtering scale (specif- 
ically, the filtering wavenumber /c_f) is defined by expanding 
the ratio of baryonic to total density fluctuations to first 
order in fc 2 : 



<5tot 



TT + TLSS 



(9) 



where fc is the wavenumber, and <5b and <5tot are the baryonic 
and total (i.e., including both baryons and dark matter) den- 
sity fluctuations, respectively. The parameter tlss (a neg- 
ative quantity) describe s the relative difference between Sb 
and <5 t ot on large scales (|Naoz fc Barkanall2007l ). i.e., 



rLSS 



<5to 



(10) 



where A = St, — <5 to t, (see also lBarkana fc Loebl l2005). The 
ratio tlss is independent of fc, and its magnitude decreases 
with time approximately oc 1/a, since A is roughly constant 
and <5tot is dom inated by the growing mode oc a (see figure 
1 top panel in lNaoz fc BarkanalbOOTT ). 

The filtering mass is defined from Hf simply as: 



M F = 



4tt _ / 1 2tt \ 3 

T P0 Ufc7J 



(ii) 



where po is the mean matter de nsity today. T his relation 
is one eighth of the definition in IGnedin! (2000) (who also 
used a non-standard definition of the Jeans mass). In fig- 
ure^] (bottom panel) we show the filtering mass (solid curve) 
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result ing from eq. (|11[) . as calculated in iNaoz fe Barkanal 
(|2007T ) (see also their figure 3). 

For each of the models we calculate the filtering mass 
as described here, assuming the model's initial conditions. 
Since the simulation is limited in box size, all of the per- 
turbations on large scales are effectively frozen in the sim- 
ulation. Therefore, we do not extract tlss directly from 
the simulations, but instead calculate it based on the ini- 
tial conditions as tlss = Ai„/(<5tot,ina), where the subscript 
"in" refers to initial. Thus, for example, for the E-8 case, 
?~lss = 0. Figure [J] (bottom panel) shows the analytical re- 
sults of the filtering mass for the fid calculation, the mean 
c s approximation and E-<5 (solid, dashed and dotted curves, 
respectively). Since the fid calculation is the most consistent 
calculation, we compare the two other models to it. 

The filtering mass represents the competition between 
gravity and pressure, as it measures the largest scale at 
which pressure has had a significant overall effect on halo 
formation. Since it measures an integrated effect over the 
formation, this mass scale is also very sensitive to the 
evolution history and the initial conditions (as shown in 
INaoz fc Barkanal 120071 ). In the mean c s model, the temper- 
ature fluct uations are greatly ove restimated on all relevant 
scales (see Naoz & Barkana 2005), while in reality the cou- 
pling to the CMB (in the fid model) keeps the temperature 
fluctuations highly suppressed for some time after recombi- 
nation. Moreover, as mentioned in section |3~T1 (and see also 
Appendix[C]), we do not include explicitly the effect of initial 
temperature fluctuations in the simulations. However, the 
temperature fluctuations from higher redshifts influence the 
baryon density at the initial redshift (see figure [TJ and sup- 
pr ess the baryon density on small scales. As demonstrated 
in lNaoz fc Barkanal l|2007h the system remembers the initial 
conditions. In other words, the initially enhanced filtering 
mass (compared to the fid model) helps maintain a higher 
filtering mass even at moderately low redshift. 

In the E-S model, the baryon perturbations start out 
much higher than in the other models, so one might expect 
that the final baryon fraction in halos would tend to be 
higher as well; here, however, it is important to separate 
two issues. The high initial baryon perturbations in the E-8 
model are present at all scales, so they affect even high-mass 
halos that are unaffected by pressure. This can explain why 
the simulation with the E-8 ICs produced the highest baryon 
fraction in high- mass halos (see the top panel of Figure [4]). 
However, when we consider the differences between large 
and small scales, the high baryon perturbations produce a 
large pressure term, increasing the effect of pressure relative 
to gravity and producing a higher filtering mass in the E-S 
model than in the fid model. Note that the filtering mass is 
particularly sensitive to the importance of pressure at the 
very highest redshifts (above 100), since at lower redshifts 
the gas cools and the Jeans mass decreases, reducing the 
contribution of these redshifts to the final filtering mass. 

We note that in INaoz fc Barkanal (2007) the calculation 
of the filtering mass in the fiducial model was compared to 
the time integrated filtering mass in a model that assumes 
the mean speed of sound model, neglects the tlss factor, 
and starts out with initial conditions as in the E-i5 model. 
Here, we have separated our discussion into several different 
cases. 



4-2.2 The non-linear characteristic mass 

There is no apriori reason to think that the filtering mass 
can also accurately describe properti e s of h ighly non-linear, 
virialized objects. For halos. iGnedinl j20 00j) defined a char- 
acteristic mass Mc for which a halo contains half the mean 
cosmic baryon fraction /;,. In his simulation he found the 
mean gas fraction in halos of a given total mass M, and 
fitted the simulation results to the following formula: 



g,calc 



= /b,0 



1+2' 



,/3 



>(£)" 



-3/o 



(12) 



where /b,o is the gas fraction in the high-mass limit. In this 
function, a higher a causes a sharper transition between 
the high-mass (constant /„) limit and t he low-mass limit 
(assumed to be f g oc M 3 ). lGnedinl (|200Ch found a good fit 
for a = 1, with a characteristic mass that in fact equaled 
the filterin g mass by his definition. By our definition, the 
clai m fromlCnedinl (12000TI is tha t M c = 8xM F . 

INaoz. Barkana fc Mesingerl ([2009) found that, given 
their errors, the filtering mass from linear theory is con- 
sistent with the characteristic mass fitted from the simula- 
tions, for two (pre-reionization) scenarios that they tested: 
the NoUV case (i.e., no stellar heating) and the Flash case 
(i.e., after a sudden flash of stellar heating). For clarity, we 
emphasize that this statement (M c = Mf) refers to our 
definition of Mf in equation (|lip . 

The characteristic mass is essentially a non-linear ver- 
sion of the filtering mass, and so it also measures the compe- 
tition between gravity and pressure. At high masses, where 
pressure is unimportant, f g — > ft,o, while the low mass tail 
is determined by the suppression of gas accretion caused by 
high baryonic pressure. 



4-2.3 Comparison between the simulation and the 
theoretical predictions 

A major conclusion of the simulation results is that differ- 
ent ICs result in different gas fractions in the final halos. 
Specifically, we measure these differences through the char- 
acteristic mass at various redshifts. varies for different ICs. 
We determine for each simulation output the characteris- 
tic mass and the parameter a using a two-dimensional fit to 
equation (|12|l . with /^o separately fixed to equal the average 
of the highest few mass bins (see Appendix|B]for a complete 
description of the fitting process, together with the 1 — a 
errors). In figure|4]we show /(,,o, a, and M c , for all the simu- 
lated cases. The characteristic mass clearly depends on the 
initial conditions, with the mean c s model and E-5 model 
both yielding gas suppression at systematically higher halo 
masses then for the fid model. The parameter a shows a less 
clear pattern with redshift, but it is generally lowest for the 
fid model. Overall, the most important implication is that 
the gas fraction in halos is highly sensitive to the assumed 
initial conditions. 

Comparing to linear theory allows us to understand 
some of these results. As noted in section 14.2.11 we calcu- 
lated the filtering mass from linear theory for each of the 
ICs, and the linear calculation allows us to understand the 
relative importance of pressure in the various IC models, at 
least during the linear evolution. Although the simulation 
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Figure 4. The parameters of the best fits in the form of equa- 
tion [TSJ different panels show M c , a, and /j, q. We consider the 
fiducial calculation, mean c s approximation and the E-<5 model 
(boxes, triangles and circles, respectively), where we fit equa- 
tion il I to all data points from halos with at least 500 particles. 
In the bottom panel we also s how the analytical calculation fol- 
lowing iNaoz fc Barkan for all the models, assuming the 
same ICs as in the simulations (solid, dashed and dotted curves 
for fid, mean c s and E-<5, respectively). We note that at 1 + z = 13 
the mean c s and the E-<5 models have the same value of M c , and 
that the fid model and the E-<5 overlap at 1 + 2 = 17. We also note 
that the data for 1 + z = 14 was unavailable due to a computer 
failure. 



seen from comparison to linear theory, the system retains the 
memory of the pressure, due to the time integrated nature of 
the filtering mass. Therefore, the higher pressure translates 
to a higher filtering/characteristic mass. 

The mean c s approximation starts with effectively 
smoother ICs than in the fid model (~ 20% underestimate of 
the small-scale baryon overdensity) . Thus, the baryonic com- 
ponents lag behind the dark matter collapse, and the pres- 
sure is always overestimated for a given baryon overdensity 
(due to the overestimated temperature fluctuations), result- 
ing in a lower gas fraction for any given halo mass, 
characteristic mass is higher than in the fid model 
cally, at z = 20 we find M c ~ 7 X 10 4 M Q and a <- 



can be compared with M c 
the fid ICs. 



3 x 10 4 M and a 



i.e., the 
Specifi- 
1. This 
0.6 for 



Recently, iHoeft et al l (|2006l ) and 

lOkamoto. Gao. fc TheunsT i 20081 ) showed th at the charac- 
teristi c mass scale does not agree with the iGnedin "fc Hui 
( 1998) filtering mass in the low-redshift, post-reionization 
regime. However, it is important to note that at these low 
redshifts, the heating/cooling and other feedback mecha- 
nisms are complex and highly inhomogeneous, so that the 
"filtering mass" calculated from linear theory is not really 
precisely defined, and the comparison of the linear and 
non-linear results canno t really be considered a direct and 
precise test. In contrast. INaoz. Barkana fc Mesingerl (|2009| ) 
found that the filtering mass gives a good approximation 
to the characteristic mas s, even in the presence of a "flash" 
heating event (see also iMesinger. Bryan fc Haimanl I2006T ) 
that is physically somewhat contrived but allows for a clear 
comparison of the linear and non-linear results. 



results come from non- linear, viralized halos, we find an ap- 
proximate agreement (typically to within ~ 20%) between 
the filtering mass, as defined here and in our previous wor k 
|Naoz. Barkana fc Mesingerl [20091 ; INaoz fc Barkanal [2007b . 
and the characteristic mass as measured in the simulation, 
for all the models. In particular, the relative sizes of M c 
among the various models, and the slow decline of all the 
characteristic masses with time, are well matched by the 
corresponding M f values predicted from linear theory. This 
close match can be understood from the fact that while both 
gravity and pressure increase during the non-linear evolu- 
tion, their relative strength only changes by a relatively 
small factor as a halo undergoes non-linear collapse and viri- 
alization. Halos in which pressure had a large effect during 
the early, linear evolution stage, keep sufficient pressure to 
maintain the suppressed baryon content all through the fi- 
nal collapse. On the other hand, in more massive halos in 
which gravity overcame pressure early on, the baryons keep 
up with the collapse of the dark matter and the pressure 
never has a major role. 

For the E-5 alternative model, we find that the resulting 
characteristic mass is higher than the result in the fid model. 
Specifically, at z = 20 we find M c ~ 5 x 10 4 M Q and a ~ 1). 
This can be understood since setting the gas fluctuations to 
be equal to the dark matter's means that the pressure of 
the gas is higher compared to the fid model. As can also be 



Summarizing our results, we find a good agreement be- 
tween the characteristic mass and the filtering mass in all the 
models. Figure [4] shows the best fitted parameters at various 
redshifts for M c and a, and our value for /b,o, for all models 
(the 1-cr (68%) confidence regions are listed in table iBl) . It 
is important to emphasize that in this statement we are re- 
ferring to our definition in equa tion (lll|) , which is one eighth 
of the original definition which IGnedin! (|2000h claimed was 
a good fit to the characteristic mass. While we have been 
careful to sel ect halos with at least 500 particle s, based on 
the results of INaoz. Barkana fc Mesingerl (|2009l ). we do not 
have the even higher mass resolution needed to perform a 
resolution convergence test as they did. Our main conclusion 
is that at least in the redshift range z = 11 — 21 the filtering 
mass provides a fairly good estimate for the characteristic 
mass. This extends the redshift range of the agreement be- 
t ween the filtering mass and t he ch aracteristic mass found 
in lNaoz. Barkana fc Mesingerl |2009l) (z = 20-25). Another 
significant result from this agreement is that previous work 
(either analytical, semi-analytical, or using simulations) that 
used the filtering or characteristic mass without accounting 
for the correct initial conditions resulted in inaccurate re- 
sults. This is due to the significant (factor of 2-3) variation 
among the predictions of the filtering/characteristic mass in 
the various models. Since this mass scale is of prime impor- 
tance in early structure formation it is imperative to calcu- 
late it accurately. 
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5 DISCUSSION 

We have used three-dimensional hydrodynamical simula- 
tions to investigate the effect of different initial conditions 
on the gas fraction in halos in the early universe. Specifi- 
cally, we studied the minimum "gas-rich" mass defined to 
have half of the mean cosmic baryon fraction. We tested 
three different models for the initial conditions (see text for 
more details) 

(i) "fi d" ICs; this mod e l is b ased on the linear evolu- 
tion from lNaoz fc Ba rkana ( 2005), which allows the baryonic 
speed of sound to spatially vary as a result of the Compton 
scattering with the CMB. 

(ii) "E-J" ICs ; in th is model, the linear evolution from 
iNaoz fc Barkanal |2005f) is modified to match a common as- 
sumption in the literature, where the linear initial overden- 
sity of the baryons is taken to be equal to that of the DM, 
i.e., 8b = 8dm = Stot, while conserving cr& from the fid model. 

(iii) "mean c s ICs"; this model assumes that the 
baryonic speed of sou nd is spatially uniform. Although 
INaoz fc Barkarlal l|2005h showed that this assumption yields 
an inaccurate evolution of the baryon density and tempera- 
ture perturbations, it is still often used in codes that gener- 
ate initial conditions for simulations. 

For all of the tests we used a total of 2 x 768 3 particles of 
dark matter and baryons with a box size of 2 Mpc, starting 
at z = 99. 

There are two major findings from the analysis we 
present here. The first, shown throughout the paper, is the 
importance of assuming the correct initial conditions, both 
for analytical calculations and numerical simulations. Struc- 
ture formation (both in the linear and non-linear regime) 
and halo gas fractions are very sensitive to the initial con- 
ditions even at relatively low redshifts (~ 10). The second 
major finding is the apparent agreement between the filter- 
ing mass and the characteristic mass (to within ~ 20%). 
This suggests, as a broader implication, that one can use 
linear theory in order to predict the overall trend of highly 
non-linear behavior (at least in the case of determining the 
gas fraction of halos). 

The the fiducial c alculation, which was presented in 
INaoz fc Barkanal (|2005T ). follows the time evolution of the 
linear overdensities correctly. However, the other ICs pro- 
duce different results for the baryonic structure formation. 
For instance, the non-linear power spectrum (fig. [3| shows 
that the system still remembers its initial condition differ- 
ences even at redshift 15. In particular, the c s model un- 
derestimates the non-linear baryonic fluctuations by about 
10% while the E-5 model overestimates them by 40% on 
small scales. 

The mean c s approximation and the E-<5 model are of- 
ten used to set the initial conditions in simula tions, e.g., the 
CMBFAST code jSeliak fc Zaldarriagal I"l996l) assum es the 
mean c s approximation while lEisenstein fc Hul (| 19991 ) is used 
with the E-<5 assumption. We have shown that the non-linear 
evolution is very sensitive to the initial conditions (figure [3]) 
and they affect the gas fraction in small halos down to red- 
shift ~ 10 (figure [4]). Our results emphasize the importance 
of the differences between the dark matter and baryons and 
of the spatial sound speed fluctuations, in both the linear 
calculation and the initial conditions of the simulations. 



It is important to emphasize that although Compton 
heating is not included in the Gadget code that we used 
in this analysis (Gadget-2), the fiducial calculation still de- 
scribes fairly well the non-linear behavior. Actually, the 
Compton heating contribution to the heating of the gas 
in non-linear objects is negligible compare to the adiabatic 
heating due to the gravitational collapse (see Appendix [Cj. 
Also, as noted above, much of the contribution to the filter- 
ing mass comes from the highest redshifts, above our simu- 
lation starting redshift of 99, since the Jeans mass is highest 
then and so the pressure has the greatest impact at that 
early time. 

In each simulation, we calculated the characteristic 
mass for which a halo keeps most of its baryons (eq. I12[) . 
We found that the fid calculation gives the lowest value, 
which suggests that with these correct ICs, pressure plays 
only a moderate role in galaxy formation. In particular, the 
characteristic mass of ~ 3 x 10 4 M is significantly below 
the minimum mass for molecular hydrogen cooling, so the 
gas content is not strongly suppressed even in the smallest 
star-forming halos. In other words we find that before sig- 
nificant heating took place the baryon fraction in halos is 
(eq. [12] with M c ~ 3 x 10 4 M Q and a ~ 0.64) 



M b ~ Mt t/b,o 



1 + 0.16 



3 x 10 4 M Q 
M 



16/25 



-75/16 



(13) 



The other alternative models give incorrect higher value 
for the characteristic mass, closer to the minimum mass 
for forming stars. Even with the fid ICs, pressure does 
strongly limit the amount of gas in minihalos below the 
molecular hydrogen cooling mass. We note that this value of 
3 x 10 4 M assumes adiabatic evolution, in particular with 
n o stellar heating. This value is con sistent with the results 
of INaoz. Barkana fc Mesingerl l|2009h for a somewhat higher 
redshift range. 

We find that the theoretical linear filtering mass (as de- 
fined in section 14.2. ip is in fairly good agreement with the 
characteristic mass. This finding is true for all the models 
tested here, throughout a significant redshift range, so this 
may imply more generally a close relation between linear 
theory and non-linear halo for mation. In addition, this is 
consis tent with the findings by INaoz. Barkana fc Mesingerl 
(2009) from AMR simulations, where the filtering mass and 
the characteristic mass agreed in the "E-<5" model, even 
when a sudden heating was introduced. 

Finally, we emphasize that our results are valid only 
in the pre-reionization era . At the end of the reionization, 
iMesinger fc Diikstral l|2008l ) concluded that the characteris- 
tic mass is likely to be close to the atomic-cooling thresh- 
old of ~ 10 s Mp), w hich is also close to the values found by 
iHoeft et al.1 (120061) andlokamoto. Gao. fc Theunsl i|200Sl ). 

Recently iTseliakhovich fc Hiratal l|2010l ) argued that 
the initial velocity difference between the baryons and dark 
matter after recombination has not been fully accounted for, 
because of a higher-order contribution that is not included 
in the linear theory approach. They estimated this higher- 
order effect within the mean c s approximation and found 
that it causes an additional suppression of the small-scale 
power spectrum, in turn affecting the formation of the first 
structures. This effect should be further investigated as in 
our detailed approach here, although this would be more dif- 
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ficult (analytically, it is a higher-order and anisotropic term, 
and to simulate it directly would require starting at quite 
high redshifts). 
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APPENDIX A: a 8 CONSERVATION 

We have defined the two different models such that they 
conserve cr%(z = 0). From linear theory we do not expect 
the evolution of the mean c s model to be significantly differ- 
ent from that of the fid model (in terms of halo abundance 
and total power spectrum). This is indeed the case for the 
evolution in time of the total fluctuations of the mean c s 
model compared to the fid model on large scales (small k), 
as shown in figure |A"T1 (lower set of thin curves). 

A more delicate treatment is needed for the E-<5 model 
(see section [2~2)l . In this case, at high redshift (such as the 
initial z — 99) , the baryons are in the process of falling into 
the DM potential. This results in a faster growth of the to- 
tal fluctuations compared to the case in which there is a 
relative velocity between the DM and the baryons (such as 
in the case of the mean cs and fid models, where the rela- 
tive velocity for the E-<5 model are negligible); see figure IaTI 
dotted thick curve. At later times, the baryon fluctuations 



approach the dark matter fluctuations, and the large scale 
behavior (i.e., on linear scales) deviates from the fid model 
by less then 0.7% (see the solid curve in fig. IA1|) . 

We also note that we have checked the overall effect 
of as(z — 0) on the main results. We have performed two 
additional simulations for the E-<5 model, where we increased 
or decreased as by 5%. We found that the calculated M c is 
within the fit errors (see appendix[B]and table lBl] ) at z > 12. 
At z ^ 12, the difference in the best fitted value is below 
0.5%." 



APPENDIX B: FIT PROPERTIES 

For each redshift snapshot for each run we find the charac- 
teristic mass and a using a two dimensional fit. In figure IBTI 
we consider two example redshifts (high, z — 19 and low, 
z = 12) for which we show the binned data points and the 
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Figure Bl. Two redshift examples of fitting the characteristic 
mass [z = 19 and z = 12). We consider the fiducial calculation, 
mean c s approximation and the E-<5 model (boxes, triangles and 
circles, respectively), where we fit equation 1121 to all data points 
from halos with at least 500 particles. We also show the fits from 
table |BT1 (dotted curves). 



resulting fit. In table IB 1 1 we show our best fit parameters. 
We note that we have checked that the fits give consistent 
results if we lower the condition on the minimum number of 
particles per halo to 300 (instead of 500). We also note that 
our determination of M c relies on an extrapolation (via the 
fit) below our simulations' resolution limit 

The parameter fb,o in equation (|12|) is an average of the 
gas fraction values in the few highest mass bins. In our simu- 
lation the high-end tail of the masses has large scatter in the 
estimated gas fraction because of the low number of halos 
(each bin among the last 3 or 4 in figure IB1I represents just 
1 or 2 halos), thus we have to average over this scatter to 
get a reasonable result. This scatter is in part a result of as- 
suming that the halos are spherical, and thus halos that are 
undergoing a major merger deviate greatly from a spherical 
shape and are treated inaccurately in our analysis. We have 
tested the resulting /^o when taking a linking parameter of 
0.1, which indeed resulted in more high-mass halos, but in 
any case was consistent with the value of fb,o we found with 
the 0.2 linking parameter. Thus, in this paper, we use the 
standard value of 0.2. 

As expected at high redshift, where we have fewer ha- 
los, the errors become quit e large . We also tried, following 
iNaoz. Barkana fc Mesingerl ([2009), to bin the data and to 
perform the fit for the binned data with the 1 — a weight 
for each bin. For the redshifts for which we had more than 
~ 1000 halos we got that the binned analysis gave results 
within the non-binned fit errors, and with comparable errors. 

We also tried the approach of taking /b,o to be a free 



parameter, but this produced very problematic fitfl This is 
mainly because of the large scatter at the high mass end, so 
that a three-parameter fit could not strongly constrain the 
parameter values. We also note the fact that /(,,o is lower 
than the mean cosmic fraction ft, by about 20% - 12% for 
the fid and mean c s models, and ~ 5% for the E-<5 model (see 
figure|4]top panel). The result in the fid and mean c s models 
may reflect the real suppression of the large-scale baryon 
fluctuations in these models; the difference i n linear theory 
is ~ 6% at z = 20 (|Naoz fc Barkana|[2007l ). but the non- 
linear evolution may increase this effect. The discrepancy in 
the E-<5 mo del may reflect a limitation of the s imulation; we 
note that in lNaoz. Barkana fc Mesingerl l|2009l ) /i,,o was also 
lower than /(, and even lower by 10% from our results at the 
overlapping redshifts (where we compare the E-5 model in 
both cases). This might be due to the fact that gas shocks 
in AMR are sharper than in Gadget simulations, and thus 
AMR may produce a more realistic gas profile, although 
the result is still below the universal cosmic baryon fraction 
(|Lin et alj|2006l ). In our simulation, going to a larger radii 
can result in a more realistic value, but we used ifeoo for 
consistency with the common definition. 



APPENDIX C: 
HALOS 



HEATING OF NON-LINEAR 



The fiducial model follows correctly the baryon density and 
temperature perturbations due to Compton scattering on 
the residual free electrons after recombination. While this 
is fully incorporated in our fid ICs, our simulation does not 
take into account Compton heating. Below we show that 
for non-linear objects the heating is actually negligible com- 
pared to the adiabatic heating due to the gravitational col- 
lapse of baryons into the dark matter potential wells. There- 
fore, it is sufficient to include Compton heating in the linear 
stage only. 

The heatin g of the gas Q C omp due to Compton heating 
from the CMB ^Naoz fc Barkanail2005T ) during the free-fall 
time \j\JGp of gravitational collapse is 



oc 4 k B (T~, 



T)pjX e (t)- 



(CI) 



where ot is the Thomson scattering cross section, p- 1 is the 
photon energy density, T 7 and T are the CMB and gas tem- 
perature and x e (t) is the electron fraction out of the total 
number density of gas particles at time t. 

The virial theorem gives a relation in collapsed objects 
between the thermal energy Eth and the gravitational energy 



E 3 
radius r„;. 



-Eg r /2 . Thus, for a halo mass M with virial 



the thermal energy can be expressed as: 



1 GM 2 



(C2) 



For all relevant redshifts and mass scales we find that 
Qcomp/Eth « 1. Therefore, neglecting the contribution of 
the Compton heating during the non-linear evolution is jus- 
tified. However, as we have shown, neglecting the Compton 



2 INaoz, Barkana fc M esingcr ( 2003) also found that treating q 
as a free parameter was unproductive. 
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heating in the linear evolution and in the initial conditions 
leads to inaccurate values for the gas fraction in halos. 



